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SUMMARY 


In  this  report  spectral  estimation,  using  Poisson  distributed  time  samples  which  have 
been  time  quantized,  is  considered.  An  expression  for  the  probability  density  function  for 
the  time  quantized  Poisson  samples  is  derived.  For  a large  class  of  nonbandlimited  spectral 
density  functions  the  bias  of  the  spectral  estimates  is  determined  and  an  asymptotic  bound 
on  the  variance  is  also  derived.  Numerical  analysis  of  the  bias  and  variance  bound  is  per- 
formed and  graphs  of  the  bias  as  a function  of  frequency  are  presented.  Finally,  examples 
of  spectral  estimates  are  given  with  varying  degrees  of  time  quantization.  Results  show 
small  amount  of  bias  introduced  into  estimates  and  variance  remains  almost  unaffected. 

The  chief  advantage  of  Poisson  sampling  theory  is  that  it  eliminates  spectral  aliasing  on 
nonbandlimited  and  undersampled  time  series  data.  Time  quantized  Poisson  data  results 
from  sampling  constraints  due  to  hardware  implementation. 
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1.  INTRODUCTION 


Spectral  estimation  is  a widely  used  method  in  the  analysis  of  ocean  data.  Typically, 
continuous  time-data  are  sampled  with  equally  spaced  data  points  and  the  spectral  density 
function  is  then  estimated  via  the  discrete  Fourier  transform  (DFT).  A new  approach  has 
been  proposed  for  certain  applications  which  use  sample  points  that  are  Poisson  distributed 
in  time  [1,2].  The  chief  advantage  of  this  approach  is  that  it  has  been  shown  to  be  alias 
free  both  by  theoretical  means  (3]  and  by  simulation  [2]. 

In  this  report  we  consider  Poisson  distributed  samples  of  a continuous  time  signal 
which  have  been  time  quantized.  Time  quantization  for  the  Poisson  case  can  be  explained  as 
follows.  The  random  nature  of  the  Poisson  distributed  sample  points  provides  a probability 
that  two  data  points  may  be  spaced  arbitrarily  close  in  time.  If  one  were  attempting  to  sam- 
ple real  world  data  with  Poisson  distributed  data  points  the  sampling  device  would  most  likely 
be  limited  in  the  speed  it  could  extract  data  samples.  Thus,  whenever  the  Poisson  distribution 
required  two  events  or  sample  points  to  come  some  small  epsilon  (e)  apart,  the  sampling  device 
would  be  unable  to  respond  in  time  to  collect  both  samples.  This  epsilon  will  be  referred  to 
as  the  “quantization.” 

There  are  several  cases  of  time  quantization  which  could  be  developed,  some  of 
which  are  listed  below. 

1 . If  an  event  occurs  during  the  quantization  time  interval  of  the  processor,  it  has 
no  effect  whatsoever. 

2.  If  an  event  occurs  during  the  quantization  time,  it  starts  quantization  time  over 

again. 

3.  If  an  event  occurs  during  the  quantization  time,  it  starts  the  quantization  time 
over  again  with  probability  p. 

4.  The  value  of  the  quantization  time  (e)  is  added  to  the  interval  after  each  event 
whether  the  event  occurs  during  the  quantization  time  or  not. 

The  first  three  cases  are  very  difficult  problems  to  investigate  analytically  because 
the  time  intervals  between  events  would  not  be  independent.  The  fourth  case  is  the  simplest 
and  is  the  focus  of  the  work  reported  here. 

In  the  following  section  the  probability  density  function  for  the  sampling  times  is 
derived.  In  addition,  using  the  spectral  estimator  proposed  in  [ 1 ] and  assuming  time  quan- 
tized sample  times,  an  exact  expression  for  the  bias  of  the  estimator  is  derived  for  a large 
class  of  spectral  density  functions.  Finally,  for  the  same  class  of  spectral  density  functions, 
an  upper  bound  on  the  variance  of  the  spectral  estimator  of  [ 1 ] is  derived  and  its  properties 
examined. 

Section  III  demonstrates  the  behavior  of  time  quantized  data  in  performing  spectral 
estimation  with  examples  of  estimation  of  a first  order  Gauss-Markov  spectral  density  using 
data  with  various  degrees  of  time  quantization.  In  addition,  the  behavior  of  the  bias  expres- 
sion is  examined  numerically  with  graphs  describing  bias  as  a function  of  frequency  for  sev- 
eral different  situations.  Finally,  in  this  section  the  behavior  of  the  variance  bound  is  exam- 
ined. Conclusions  and  a discussion  of  the  time  quantization  problem  are  included  in 
section  IV. 
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II.  THEORY 


In  this  section  a derivation  of  the  probability  density  function  for  the  time  quantized 
sampling  times  is  presented.  Also  determined  are  the  bias  and  a bound  for  the  variance  of  the 
spectral  estimator  of  1 1 ) for  time  quantized  Poisson  sampled  data. 

Let  X = (X(t),  - °°<t<°°}  be  a continuous  parameter  zero  mean,  second  order, 
stationary  random  process  with  correlation  function  C(t)  e L|  and  spectral  density  function 
ip(A)  given  by 


= 


C(t)e 


-iAtil 

2ir 


(1) 


The  time  quantized  sampling  instants  are  driven  by  a Poisson  point  process  and  are 
given  by 

t0=  0 

tn  = tn-l+an  + e n = 1,2,... 

where  the  {an } are  independent  identically  distributed  random  variables  with  a common 
exponential  distribution  F(x)  = 1 - e“0x.  It  is  assumed  that  the  sampling  instants  {tn  }“_q 
are  independent  of  the  process  X.  Note  that  (3  is  the  average  sampling  rate  and  as  shown 
in  1 1 1 , almost  every  realization  of  a Poisson  point  process  has  a density  of  points  0.  The 
effect  of  the  quantization  (e)  on  the  sequence  {tn  } is  to  enlarge  each  interval  by  e and  there- 
fore shift  the  n1*1  sample  time  tn  by  ne. 

Given  the  observation  {X(tk)}k=|  we  will  consider  the  estimate  ^(A)  of  the  spec- 
tral density  <p(A)  as  proposed  in  [ 1 1 : 

MN  N-n 

*N(X,“Sn2  £ X(tk+n)X(tk)cosA(tk+n-tk).  (2) 

n=  1 k = l 

where  Mn  > 0 is  such  that  Mn  00  and  Mn/N  0 as  N -*• 

The  following  Lemma  will  be  needed  in  the  proof  of  the  bias  and  variance  of  the 
time  quantized  spectral  estimate  (2). 


Lemma  1:  Consider  the  probability  density  function  for  the  distribution  of  the  intervals  of 
a Poisson  process. 

fs(s)  = 0e~0s  U(s) 


where  U(  * ) is  the  unit  step  function.  If  the  probability  density  function  for  intervals  of  a 
Poisson  process  that  have  been  quantized  or  time  expanded  by  an  amount  e is 

fs(s)  = 0e-0(s-e)U(s-e)  , 


then  the  probability  density  function  for  the  time  quantized  Poisson  distributed  sampling 
times  {tn } is  given  by 


ftn(t)=0n 


(t  - en)n_* 


e-0(t  - en)  t _ en) 


(n-  1)! 


(3) 


Proof:  Since  the  intervals  {si } of  a Poisson  process  are  independent,  the  sampling  times{tn} 
can  be  represented  by  a sum  of  independent,  identically  distributed  random  variables.  If 
each  interval  is  quantized  by  an  epsilon,  we  then  have 


t2  = si  +s2+2e 
t3  = si  +s2  + s3  + 3e 


n 


ln  = X S> 


t„  = > sj  + ne  . 
i=l 


Now  the  characteristic  function  of  the  time  quantized  intervals  can  be  expressed  as: 
<£t(co)  = / U(s  - e)  P e-0(s-e)  e'^8  ds 

•Loo 

= (^—)ek“ 

\0  - 1 ojJ 

The  characteristic  function  of  the  n1*1  sample  time  will  then  be: 

<Ptn(u)  = [ «^t( oo) ] n = eiewn. 

Finally,  the  density  function  ftn(t)  is  the  inverse  Fourier  transform  of  ^tn(co) 

r oo 

/-oo  ' ' 

= gri  (t~enf"~1  e~^(t  ~cn)  u(t  - en).  e>0 

(n  - 1)! 


Q.E.D. 

It  is  easily  seen  that  ftn(t)  satisfies  the  properties  of  a probability  density  function 
1 4 1 namely 

ftn(t)  > 0 V t 
/ ftn(t)dt=l 

Z— OO 
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Figure  1 is  a graph  representing  the  probability  density  function  ftn(t)  as  a function 
of  time  for  several  values  of  n. 


STATISTICS  OF  THE  SPECTRAL  ESTIMATOR  FOR  THE  TIME  QUANTIZED  CASE 

Here  we  assume  that  the  process  X is  a stationary,  Gaussian  process.  For  a large  class 
of  spectral  densities  we  will  present  an  expression  for  the  bias  and  a bound  for  the  variance 
of  the  spectral  estimator  (2)  for  time  quantized  sample  points. 

Let  the  covariance  function  C(t)  have  the  following  form 

C(t)  = A e-“|T|  (4) 

with  corresponding  spectral  density  function 


where  A and  a are  Fixed  positive  constant  parameters,  and  define  the  bias  of  the  estimate  (2) 
to  be 


b[£N(A)]  = EI^IA)]  -tp(\). 


Theorem  1 : Let  the  covariance  function  C(t)  be  given  by  (4),  then  the  bias  of  the  estimate 
(2)  for  the  time  quantization  e is  equal  to: 


mn 

b[^N(X)1  ^ (l  -£j(0e^a)n  [(a  + 0)2  + X2fn/2 


• cos 


(Aen)  + n tan-* 


(»*«)  ‘'(,2^2)' 


(6) 


Proof:  Let  Ex  and  E { tn } represent  expectation  with  respect  to  X and  { tn } . Then 


b[*N(A)]  = E{tnJEx 
M 


M 


£ X(tk+n>X(tk)cosX<tk+n-tk) 


N N-n 


n=l  k=l 


-^A) 


A 

7T0N 


N r° 

y (N-n)  / 
n=l  Jo 


C(t)  cos  (At)  ftn(t)  dt 


"*(a2  + A2) 


where  ftn(t)  is  given  by  (3).  Substituting  for  C(t)  and  ftn(t)  we  then  have 
MN  f 00  n_i 

•JjI<N-n)  / e^'“!(M)57^Tnrjne"Wt'en),“ 

n=  1 Jne 

® V2  + A2/ 
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figure  1.  Probability  density  function  for  time-quantized  Poisson  distributed  events  where  0 = 10  rad/sec,  and  e = 0.01. 


M 


= ^ (N-n) — - e^en  / e"*u  en>'^  + cos  X(u  + en)  un“*  du 

7T0N  Z-i  (n  - 1 )!  / 

n=l  J 0 

"•(tTt) 

by  a change  of  variable  u = t - ne.  Next  expand  the  cosine  function  to  give 

/*oo 

(N-n) — - -e^en  e-en((3  + a)  cos  (Xen)  / e_u^  + a!^  cos  (Xu) 

7T/3N  (n  - 1 )!  / 

n=  1 J0 


• un  ' du  - sin  (Xen)  / e 


u(^  + a)  sin  (Xu)  un  ■ du  - — (—  a — 

v \cx2  + \2) 


mn 

=— y 

ttON  ZL 


(N  - n) 


Pn 


3-en  a 


tt0N  Z^  ' (n  - 1)! 

n=l 


cos  (Xen)  F(n)  |(a  + P)2  + X“]  ^ 


7T  ft 


n tan  1 I I 

V*  + /V 

["  (^)] 

mn  , . 

20-s) 


- sin  (Xen)  F(n)  [(a  + (5)2  + X^| 


A a 


-n/2 


77  a2  + X“ 


(fie-*01)  [(or  + /3)“  + X“1 


-n/2 


_1  / A \ A/  a \ 

cos  Xen  + n tan  *1 — — ) - - ( r-l 

\a  + Pj]  77  \a“  + X“/ 


where  the  integrals  were  evaluated  using  [5]. 

q.e.d. 


An  examination  of  (6)  reveals  that  since  the  first  few  terms  in  the  summation  are 
dominant  there  is  not  a high  dependence  on  N and  Mn.  We  see  also  that  the  bias  decreases 
as  X increases  at  a rate  of  at  least  1 /X. 

Next  we  consider  the  variance  of  the  estimate  (2).  We  assume  that  the  process  X is 
stationary  to  order  four  and  has  covariance  and  spectral  density  functions  as  given  in  (4)  and 
(5),  respectively. 

In  the  derivation  of  the  variance  bound  we  begin  with  the  same  procedure  as  in  [ 1 ] . 
Consider  the  fourth  order  cumulant  function  defined  by 


0(r , , r2>  r3)  = E[  X(t)  X(t  + r , ) X(t  + r2)  X(t  + r3)) 

-C(tj)C(t3 -r2)-C(r2)C(T|  - t3)  - C(t3)  C(r2  - r j ) (7) 

where  X and  C(r)  are  as  previously  defined. 
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If  we  assume  the  process  X to  be  Gaussian,  then  (7)  becomes 

E[X(t)  X(t  + r ] ) X(t  + r2)  X(t  + 73)]  = C(t,  ) C(r3  - r2)  + C(r2)  C(t,  - 73) 

+ C(73)C(72-7j). 

Now,  we  look  at  the  second  moment  of  the  estimate  (2) 


MN  MN  N-n  N-n' 


C(tk+n  "‘k^^k'+n'  _tk') 


(Wn-ln 

+ C(‘k'  - lk>  C(tk+n  - tk'+n')  + C(tk-+n-  - tk)  C(tk'  - tk+n) 
• cos  X(tk+n  - tk)  cos  Mtk'+n'  - tk'). 


(8) 


The  expectation  with  respect  to  tn  in  (8)  makes  its  solution  very  difficult  to  obtain.  In  order 
to  solve  the  expectation  in  (8)  and  bound  the  variance  of  (2)  we  will  break  the  summations 
up  and  evaluate  them  in  eleven  disjoint  regions  corresponding  to  the  different  permutations  of 

{k,  k',  k + n,  k'  + n' }. 

Let  E [ ( X )]  ~ be  divided  into  the  following  three  expressions 
3 

E[^N(X)]2  = ^ Uj(X) 
i=l 


where 


Uj(X)  = 


MN  MN  N-n  N-n' 

2 2 2 2 


n=ln’=l  k=l  k'=l 

• cos  Mtk+n  ' V cos  ^k'+n'  “ lk’> 
MN  MN  N-n  N-n' 


C(tk+n  ~ lk^  C^lk'+n'  ~ lk'^ 


U2(X)=—L-^^  S ^ S Etn 


(77/3N) 

cos  X(tk+n  - tk)  cos  X(tk-+n»  - V) 
MN  MN  N-n  N-n' 


c(tk'  - tk)  C(tk+n  - tk-+n') 


u3(X)  = 


2 2 2 2 b. 


(^N>Z  n=l  n - 1 k=l  k*l 


^lk'+n'  ~ lk^  ^(tk'  ~ *k+n* 


• cos  Mtk+n-tk)cosX(tk»+n--tk') 


(9) 


GO) 


(11) 
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Lemma  2:  If  the  process  X is  Gaussian  and  has  a covariance  function  given  by  (4>  then  an 
upper  bound  on  Uj(a.)  is  given  by: 


2 Mn 

U,(X)  < Y,(X)  = — — — - ^ (N-n)W2n  + 


2A‘ 


(7r/3N)^ 


IMvi  Mxi  < 

4Z  Zj  \[(ce  +/3)  - i\]  / RC\j(a  + 0)-i\l/ 

j^(N-n)(N-n'-n+  1 ) - i(N  - n)(N  - n + 1 )J  I(N  - n - n'  - 1 ) 

mn  mn 

(N-n')  Wn'(n'  - 1)+  J 


n =1 


’ MN  N-n  N-n'-k 

1,1,  I ,wstn' 

n=l+n'  k=N-n'-n  s=n-n' 


M 


N n+1 


♦ II  (N-n'-n-  l)Ws+n' 

n=l+n' s=n-n' 

MN  r 2 -.1 

+ £ (N-n)Wn  |^-+2-|nJ  j 


where 


W = 


-ea 


a + 0 

Proof:  See  appendix. 

Lemma  3:  If  the  process  X is  Gaussian  and  has  a covariance  function  given  by  (4)  then  an 
asymptotic  bound  of  l^fX)  as  N -*•  °°  is  given  by: 

«*»  «=  v2,x, . > U Re  L t „ . ps±a,\l 

(np)2  N [2  |_  \(c*  + 0)2  + X2/ 

mn  mn  N_n<_, 

£ I I (N  - n'  - s)  Ws+n'  ( (N  - n - n - 1 ) 

" n=I  n'=I  s=n 


(a+0-iX) 

2[(a+0)2  + X2 1 2(a+0)cos(eX)-2Xsin(eX)] 

(<*+|3-iX)-|3e_e(a'iA) 

( (a+0)2  + X2l  -0e_eQ![  2(cr+0)cos(eX)-2Xsin(eX)l+02e“2cQ: 
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<2 


r 


2 (a  + P)  - 2ffe  eQ!  cos  (ea) 


[(a  + 0)2  + X2]  - 2j3e_CQ![(a  + 0)  cos  (eX)  - X sin  (eX)J 
0-iX 


+ 02e_2ea 


+ Rt  r £±tTarwir-2TfiLilt  « [i-  ♦ (-J^rl T \ 

^(a  + /3-iX)-)3e  e(a!  lX)Jj^(a  + /3)2  + X2J  \(a+ 0)2+X2  / 

+2R«r  “+,,-ix 


I (o  + jS-iX)-^-6^-^) 


• 2 

I 


Proof:  See  appendix. 

Lemma  4:  If  the  process  X is  Gaussian  and  has  covariance  function  (4)  then  an  upper  bound 
on  U3(X)  is  given  by: 


M 


N 


2AZ  V’ 

U3(X)<Y3(X)  = > 

(7T/3N)2  ", 


’mn 

„tA(w-D2  w-‘  ) 


M 


N N-n 

2 2 

n'=n  k=N-n'-n+2 


M 


W 


W-  1 


+ Wnj^3n2N  - 


6nN  + 6n2  - 


N-k-W')t  Y ,N-n'-nH)Wn|Wn;n~'-1> 

„r?„  wn  - i 

-t] 


+ (N  - n)  W2n  + [4N  - 2]  Wn 


2A‘ 


mn 

2 


(ir/3N)z  n»_, 


■mn 

V 1N-n--.^lW-(WnT 

^ Wn  - I 


n=n 


M 


N N-n 

2 2 


W 


W - 1 


(WN~k+l  _wn')  + (N-n')(n'  - 1)  Wn' 


n=n'  k=N-n-n'+2 
Note  that  Y3(-)  does  not  depend  on  X. 

Proof:  See  appendix. 


Theorem  2:  If  X is  a continuous  parameter  zero  mean,  fourth  order,  stationary,  and  Gaussian 
random  process  with  covariance  function  given  by  (4),  then  an  asymptotic  bound  on  the 
variance  of  the  estimator  in  (2)  is  given  by 
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(12) 


Var  [^jsj(X))  = U,(X)  + U2(X)  + U3(X)  - E^[^j(X)) 

< |Yj(X)|  + |Y2(X)|  + |Y3(X)|  - |E2l^N(X)]  1 


as  N -*•  °°,  Mjyj  -*  °°,  Mjyj/N  -*■  0. 

Proof:  The  proof  follows  from  Lemmata  2,  3,  and  4. 

The  forms  of  YjCX),  Y2(X)  and  Y3(X)  are  unfortunately  very  awkward  and  as  a 
result,  it  is  not  easy  to  see  how  these  equations  behave.  The  summations  in  Y j(X)  and  Y3(X) 
may  be  eliminated  through  the  long  and  tedious  procedure  of  using  geometric  expansion 
type  identities.  This  is  presently  being  performed  and  when  accomplished  will  greatly  sim- 
plify the  variance  expression  and  provide  an  insight  to  the  behavior  of  the  variance.  In  addi- 
tion, if  the  bounds  Y j(X)  and  Y3(X)  can  be  shown  to  be  maximized  at  e = 0,  then  they  can 
be  bounded  by  the  results  in  [ 1 ) : 

UjCX)  - E2[^n(X)1  = 0(1/N) 

U3(X)  = 0(1 /N) 
uniformly  in  X as  N -*■  °°. 
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III.  NUMERICAL  AND  SIMULATION  RESULTS 

Since  the  theoretical  analysis,  in  Section  II,  showing  the  bias  and  variance  of  the  time 
quantized  spectral  estimator,  is  somewhat  involved  and  not  easily  interpreted,  this  section 
will  provide  insight  into  their  behavior.  Here  we  will  present  graphs  of  the  bias  as  a function 
of  frequency  for  several  e//3  ratios,  a numerical  evaluation  of  the  variance,  and  examples  of 
time-quantized  spectral  estimates  using  simulated  time  series  data. 

Figure  2 is  a graph  of  the  spectral  density  for  the  process  X where  covariance  function 
parameters  a and  A have  been  set  equal  to  1 . This  is  representative  of  a nonbandlimited, 
low-pass  spectrum.  The  peak  of  this  spectrum  is  at  X = 0,  the  half  power  point  at  X = 1 , and 
the  spectrum  decreases  at  a rate  proportional  to  l/X^.  By  adjusting  the  parameters  a and  A, 
the  shape  of  the  spectrum  can  be  changed  to  broaden  or  narrow  the  peak.  This  spectral 
density  will  be  used  in  the  following  numerical  analysis  of  the  bias  and  variance.  The 
basic  form  of  the  covariance  function  C(t)  used  here  which  results  from  the  solution  of  a 
first  order  homogeneous  differential  system  can  be  combined  [6]  by  addition  to  form  the 
solution  to  a homogeneous  differential  system  of  arbitrary  ordei. 


BIAS. 


Let  the  bias  of  the  spectral  estimator  (2)  be  given  by  (6).  We  then  define  the  per- 
cent bias  as: 


% bl^X)]  = 


b[^,(X)] 

-MX) 


Figures  3 through  8 are  graphs  of  both  the  bias  and  percent  bias  as  a function  of  fre- 
quency X for  several  values  of  the  average  sampling  rate  (3  and  time  quantization  e.  Three  average 
sampling  rates  were  used  for  the  calculations  of  the  bias.  Figures  3 and  4 were  generated  with 
0 = 1 rad/sec.  Figures  5 and  6 were  generated  with  0 = 5 rad/sec  and  Figures  7 and  8 were  generated 


Figure  2.  Spectral  density  function  simulated  for  numerical  analysis. 
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Figure  6.  Percent  bias  of  time-quantized  spectral  estimate  as  a function  of  frequency  and  epsilon  for  /l  = 5 rad/sec,  N = 1000  and  M, 


time-quantized  spectral  estimate  as  a function  of  frequency 


LEGEND 


Figure  8.  Percent  bias  of  time-quantized  spectral  estimate  as  a function  of  frequency  and  epsilon  for  0 = 10  rad/sec,  N = 1000  and  M 


with  /3  = 1 0 rad/sec.  Each  of  the  graphs  contains  4 curves.  These  curves  represent  time  quan- 
tization epsilons  of  0.0 1 , 0. 1 , 1.0  and  1 0.0  percent  of  the  average  sampling  rate  (3.  The  hori- 
zontal axis  of  each  graph  represents  frequency  and  the  bias  and  percent  bias  are  computed 
at  frequency  increments  of  1 . The  value  of  used  is  6 for  /3  = 5,  10  rad/sec  and  = 7 
for/3  = 1 rad/sec. 

The  behavior  of  the  bias  for  the  larger  values  of  epsilon  is  oscillatory  and  increases 
with  sampling  rate.  The  values  of  the  bias,  although  larger  for  X < 1 0,  are  overall  fairly  small 
but  as  expected,  increase  as  e becomes  larger.  This  is  probably  due  to  spectral  aliasing  of  high 
frequency  components. 

In  order  to  show  the  relationship  of  the  bias  to  the  spectral  density,  Figures  4,  6, 
and  8 represent  percent  bias.  The  percent  bias  is  the  percent  of  the  amplitude  of  the  spectral 
density  which  is  offset  due  to  bias.  The  parameter  values  in  these  three  graphs  are  the  same 
as  for  the  bias  graphs.  Each  of  Figures  4,  6 and  8 shows  a percent  bias  which  increases  with 
epsilon.  The  graphs  also  show  a divergence  of  percent  bias  as  frequency  increases.  This  can 
be  misleading  unless  interpreted  correctly.  It  should  be  clear  that  the  spectral  density  becomes 
smaller  as  frequency  becomes  larger  which  means  that  a constant  value  of  bias  becomes  a 
larger  percent  of  the  spectral  density  as  frequency  increases.  These  graphs  also  show  that  the 
amplitude  of  the  spectral  density  decreases  at  a rate  greater  than  the  bias. 

VARIANCE 

The  asymptotic  variance  bound  for  the  time-quantized  spectral  estimates  given  in 
( 1 2)  is  very  cumbersome.  There  is  little  that  can  be  said  about  it  in  its  present  form.  There- 
fore, in  order  to  provide  some  insight  to  its  behavior,  it  has  been  evaluated  numerically  for 
several  sets  of  parameters.  Two  values  of  epsilon  and  two  sampling  rates  were  chosen  in  the 
numerical  evaluation  with  N = 1000  and  = 6.  Table  1 shows  the  bound  on  the  variance 
of  £n(X)  along  with  V j ( A ),  Y2C  A),  V3(X)  for/3  = 1 and  e//3  = 0.0001 . The  values  of  the  inter- 
mediary functions  Y j ( A),  Y->(X),  Y^A)  are  shown  in  order  to  give  the  reader  a feeling  for  their 
relative  contributions  to  the  bound  of  the  variance  (12).  Note  that  the  major  contribution 
to  the  variance  comes  from  Y2(X)  as  predicted  in  1 1 1.  In  addition,  the  variance  is  almost 
invariant  with  X except  at  X = 0.  Tables  2,  3,  and  4 provide  numerical  results  for  /3  = 1 rad/sec, 
e = 0.01  :p  = 10  rad/sec,  e = 0.001  ;and  /3  = 1 0 rad/sec,  e = 0.1,  respectively. 


Table  1 . Numerical  results  showing  second  moment  component  bounds  of  spectral  estimator 
i^N(X)  and  the  asymptotic  variance  bound  of  estimator  for  0=  1 rad/sec,  e = 0.0001 . 


X 

Y ](A) 

Y2(X) 

Y3(X) 

Var  (£f\|(X))  bound 

0 

2.401  X 10~2 

0.243 

4.065  X 10"3 

1.733  X 10"1 

10 

5.313  X I0-4 

4 952  X 10‘2 

4.065  X 10"3 

5.411  X 10~2 

20 

5.326  X 10'4 

,45  X 10-2 

4.065  X 10"3 

5.305  X 10'2 

30 

5.331  X 10'4 

4.825  X 10"2 

4.065  X 10'3 

5.285  X 10"2 

40 

5.333  X 10'4 

4.81  8 X I0~2 

4.065  X 10"3 

5.278  X 10'2 

50 

5.333  X 10'4 

4.815  X 10"2 

4.065  X 10'3 

5.275  X 10-2 

60 

5.333  X 10‘4 

4.813  X I0~2 

4.065  X 10'3 

5.273  X 10"2 

70 

5.333  X 10~4 

4.81  2 X I0-2 

4.065  X 10‘3 

5.272  X 10'2 

80 

5.332  X 10"4 

4.81 1 X 10-2 

4.065  X I0'3 

5.271  X 10"2 

90 

5.331  X I0'4 

4.81 1 X 10-2 

4.065  X 10-3 

5.271  X I0‘2 

100 

5.330  X 1 O'4 

4.81 0 X 10‘2 

4.065  X 10'3 

5.270  X 10'2 

20 


Table  2.  Numerical  results  showing  second  moment  component  bounds 
of  spectral  estimator  i£n(X)  and  the  asymptotic  variance  bound 
of  estimator  for  0 = 1 rad/sec,  e = 0.0 1 . 

X 

Y j(X) 

Y2(X) 

y3(X) 

Var  (<?n(X))  bound 

, 0 

2.400  X 10~2 

0.239 

3.942  X 10-3 

1.728  X 10-1 

10 

5.313  X I0'4 

4.857  X 10"2 

3.942  X 10'3 

5.304  X 10"2 

20 

5.326  X 10'4 

4.752  X 10~2 

3.942  X 10'3 

5.199  X 10'2 

30 

5.331  X 10-4 

4.733  X 10~2 

3.942  X 10'3 

5.180X  I0'2 

40 

5.333  X 10*4 

4.727  X 10~2 

3.942  X 10"3 

5.174  X 10"2 

50 

5.333  X 10-4 

4.725  X 10~2 

3.942  X 10‘3 

5.172  X 10"2 

60 

5.333  X I0-4 

4.725  X 10~2 

3.942  X 10"3 

5.172X  10'2 

70 

5.333  X 10~4 

4.726  X 10‘2 

3.942  X 10"3 

5.173  X 10'2 

, 80 

5.332  X 10'4 

4.727  X 10'2 

3.942  X I0-3 

5.174  X 10"2 

1 90 

5.331  X 10~4 

4.729  X 10-2 

3.942  X 10-3 

5.1  76  X 10"2 

100 

5.329  X 10~4 

4.732  X 10-2 

3.942  X 10"3 

5.1  79  X 10~2 

Table  3.  Numerical  results  showing  second  moment  component  bounds 
of  spectral  estimator  <£n(X)  and  the  asymptotic  variance  bound 
of  estimator  for  (3  = 10  rad/sec,  e = 0.001 . 

X 

Y ](X) 

Y2(X) 

Y3(X) 

Var  (Jn(X)) 

4.779 

X 

10' 

6.187 

X 

10' 

6.193 

X 

10' 

6.187 

X 

10' 

6.187 

X 

10' 

6.187 

X 

10' 

6.188 

X 

10' 

6.188 

X 

10' 

6.188 

X 

10' 

6.188 

X 

10' 

6.188 

X 

10' 

3.553 

X 

10' 

4.435 

X 

10' 

4.598 

X 

10 

4.654 

X 

10' 

4.680 

X 

10 

4.693 

X 

10 

4.702 

X 

10' 

4.707 

X 

10' 

4.71  1 

X 

10' 

4.714 

X 

10' 

7.626 

X 

10' 

7.626 

X 

10' 

7.626 

X 

10 

7.626 

X 

10 

7.626 

X 

10' 

7.626 

X 

10 

7.626 

X 

10' 

7.626 

X 

10' 

7.626 

X 

10' 

7.626 

X 

10' 

7.626 

X 

10 

4.506 

X 

10' 

4.467 

X 

10 

4.680 

X 

10 

4.736 

X 

10 

4.762 

X 

10 

4.775 

X 

10' 

4.784 

X 

10 

4.789 

X 

10' 

4.793 

X 

10' 

4.796 

X 

10' 

Table  4.  Numerical  results  showing  second  moment  component  bounds 
of  spectral  estim  ..or  <£n(A)  and  the  asymptotic  variance  bound 
of  estimator  for  (3  = 10  rad/sec,  e = 0.1 . 


X 

Y i (X) 

Y2(X) 

y3(X) 

Var  (£n(Xj) 

0 

2.604  X 10"3 

1.550 

3 403  X lO’4 

1.543 

10 

6.441  X 10‘5 

2.109  X 10"2 

3.403  X 10^ 

2.139  X 10~2 

20 

5.810  X 10"  5 

2.509  X 10'2 

3.403  X I0-4 

2.541  X 10'2 

30 

5.427  X 10'5 

3.344  X 10~2 

3.403  X I0*4 

3.381  X 10‘2 

40 

4.002  X 10'5 

5.624  X 10"2 

3.403  X 10^ 

5.661  X 10"2 

50 

4.942  X 10"5 

7.090  X 10"2 

3.403  X 10"4 

7.124  X 10~2 

60 

4.853  X 10'5 

5.252  X I0'2 

3.403  X I0-4 

5.291  X 10'2 

70 

3.925  X !0*5 

4.1 13  X 10"2 

3.403  X 10-4 

4.150  X 10'2 

80 

4.046  X I0"5 

3.890  X I0‘2 

3.403  X I0-4 

3.927  X I0‘2 

90 

3.867  X 10'5 

4.276  X I0"2 

3.403  X 10-4 

4.314  X 10‘2 

100 

3.820  X 10"5 

5.130  X 10"2 

3.403  X 10-4 

5.168  X 10"2 

For  the  sake  of  illustration  we  generated  time  series  data  as  in  [2]  for  the  spectral 
density  of  (5)  with  A = a = 1 , sampled  it  using  several  values  of  epsilon  and,  using  the  esti- 
mator (2),  computed  spectral  estimates.  The  results  are  shown  in  Figures  9,  10,  and  1 1 for 
P = 1,5,  and  10  rad/sec,  respectively.  Here  the  epsilons  used  are  the  same  as  those  used  in 
computing  the  bias  and  percent  bias.  In  these  graphs  the  true  spectral  density  ^>(A)  is  repre- 
sented by  a solid  curve  without  symbols.  The  vertical  axis  represents  amplitude,  while 
frequency  is  along  the  horizontal  axis.  Any  negative  estimates  generated  are  set  equal  to 
zero  since  the  spectral  density  by  definition  is  everywhere  positive.  The  number  of  data 
points  used  was  N = 1000  and  was  chosen  by  generating  spectral  estimates  using  several 
values  of  and  selecting  the  which  provided  the  best  estimate  (based  on  mean  squared 
error). 

In  Figure  9 where  the  lowest  average  sampling  rate  was  used  (P  = 1 rad/sec),  the  peak 
of  the  spectral  density  is  generally  estimated  quite  well.  As  the  average  sampling  rate  is 
increased  (Figures  10  and  1 1 ),  the  values  for  the  estimates  of  the  peak  become  biased.  This 
bias  is  predicted  in  Figures  3,  5,  and  7 where  bias  near  the  origin  increases  both  as  a function 
of  e and  as  a function  of  p.  It  also  appears  in  Figures  9,  10,  and  1 1 that  spectral  estimates 
for  X > 5 improve  with  increasing  the  value  of  /3.  There  is  no  periodic  type  structure  in  any 
of  the  estimates  that  would  suggest  spectral  aliasing  as  one  might  expect,  especially  for  the 
larger  epsilons. 

Spectral  estimates  for  the  case  when  e = 0 were  computed  and  observed  to  be  almost 
identical  for  the  case  e/P  = 0.01%.  For  this  reason,  spectral  estimates  for  e = 0 are  not 
shown  on  the  graphs. 
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a = EPSILON' 0.0001 
O = EPSILON  = 0.0010 
A ' EPSILON '0.0100 
+ = EPSILON' 0.1000 
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□ = 

EPSILON  - 

0.0010 

o * 

EPSILON  = 

0.0100 

A = 

EPSILON  = 

0.1000 

+ a 

EPSILON  ■= 

1.0000 

1 ' ■ ■ T ■ ■ T" " 


0.0  2.0  4.0  6.0  8.0  10.0  12.0  14.0 

LAMBDA 


Figure  1 1 . Spectral  estimate  for  0 = 10  rad/sec,  = 6 and  several  values  of  epsilon. 
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IV.  CONCLUSIONS 


This  report  has  shown  results  of  an  examination  of  the  spectral  estimator  proposed 
in  ( 1 ] for  time-quantized  Poisson  distributed  sampling  intervals.  The  probability  density 
function  for  the  time-quantized  Poisson  samples  as  well  as  the  bias  and  a bound  for  the 
variance  of  the  spectral  estimator  for  a large  class  of  spectra  have  been  derived. 

On  performing  a numerical  analysis,  the  bias  of  the  spectral  estimator  was  seen  to 
increase  with  increasing  average  sample  rate  (3  and  decrease  with  frequency  X at  a rate  be- 
tween X and  X“.  In  all  cases  examined  the  bias  appeared  to  be  largest  in  the  region  around 
X = 0 and  increased  with  epsilon. 

The  bound  for  the  variance  of  the  spectral  estimator,  which  was  derived  through  a 
complicated  procedure,  resulted  in  a cumbersome  asymptotic  expression.  With  useful 
information  difficult  to  extract  from  this  bound,  it  was  evaluated  numerically  and  found 
to  be  small  with  respect  to  the  spectral  density  except  near  the  origin.  In  addition,  except 
near  the  origin,  it  was  found  to  have  small  dependence  on  the  value  of  the  average  sampling 
rate. 

It  was  shown,  using  synthetically-generated  data,  that  the  behavior  of  the  time 
quantized  spectral  estimator  was  similar  to  that  predicted  by  the  theoretical  bias  particularly 
near  the  origin  and  not  adverse  to  that  of  the  variance  bound.  This  simulation  also  showed 
that  degradation  of  the  spectral  estimates  was  small  for  e/0  ratios  of  less  than  0.1%,  but  for 
larger  e//3  ratios  the  effect  was  quite  pronounced.  The  measure  of  degradation  used  was 
integrated  mean  squared  error  over  the  range  of  the  graphs.  It  is  not  clear  if  mean  squared 
error  is  the  best  measure  of  degradation  since  other  measures  were  available,  but  it  was  used 
in  order  to  be  consistent  with  |2).  Another  measure  of  degradation  which  could  be  used  is 
the  quality  of  estimates  of  the  peak  of  the  spectrum.  For  this  measure  the  results  are  the 
same;  larger  epsilon  yields  more  degradation. 

The  results  in  this  report  indicate  that  for  the  class  of  spectral  densities  examined,  a 
time  quantization  of  the  Poisson  sample  point  of  up  to  0.1%  of  the  average  sampling  rate 
may  perhaps  be  tolerated. 

There  are  several  questions  which  remain  unanswered  concerning  this  time  quantized 
Poisson  sampling  problem  and  should  be  pursued.  First  of  all,  the  asymptotic  bound  for  the 
variance  of  the  estimator  should  be  simplified  and  a determination  of  its  tightness  be  made. 

It  is  not  clear  that  for  e = 0 the  bias  and  variance  results  reported  here  are  bounded  by  or 
equal  to  the  results  of  [ 1 ].  This  should  be  determined.  The  effects  of  using  different 
covariance  weighting  functions  should  be  examined  for  the  time  quantized  case.  Finally, 
other  more  complicated  spectral  densities,  such  as  a narrowband  signal  modulated  by  a 
carrier  signal,  could  be  examined  for  results  comparable  to  those  reported  here. 
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V.  APPENDIX 


Here  we  present  proofs  of  Lemmata  2,  3,  and  4 which,  when  combined,  make  up  a 
bound  for  the  second  moment  of  the  spectral  estimator  (2)  for  the  Gaussian  case,  with  zero  mean. 

As  was  previously  mentioned  in  section  II,  the  expectation  with  respect  to  tn  in  (8) 
cannot  be  evaluated  unless  it  is  examined  over  regions  involving  the  different  permutations 
of  {k,  k',  k + n,  k'  + n' }.  These  nonintersecting  regions  Rj  are  defined  in  [ 1 ] and  are: 

Rj  =(k<k  + n<k'<k'  + n' |k,k'6R) 

R:  = (k  < k'  < k + n < k'  + n'  1 k,k'  € R) 

R3  = (k  = k'  < k + n < k'  + n'  | k,k'  € R) 

R4  = (k  < k'  < k + n = k'  + n'  | k,k'  € R) 

Rc  =(k  = k'<k  + n = k'  + n' |k,k' €R) 

(A-l) 

R6  = (k  < k < k + n'  < k + n | k,k'  £R) 

R?  =(k  = k'<k'  + n'<k  + n |k,k'  € R) 

R8  =(k'<k'  + n'<k<k  + n | k,k' €R) 

R9  =(k'<k<k'  + n'<k  + n | k.k' €R) 

R,0=(k'<k<k'+n'  = k + n | k,k'  €R) 

R , , =(k'<k<k  + n<k'  + n'  | k,k'  £ R) 

The  union  of  these  regions  is  equivalent  to  summing  over  the  indices  k and  k'  of  (8) 
or 

1 1 

R = U Rj. 
i=l 

There  is  considerable  symmetry  between  pairs  of  the  regions  defined  in  (A-l ) which 
can  be  used  in  evaluating  (9),  (10),  and  (11).  Let  Uk(X),  k = 1 , 2,  3 be  the  sum  of  Uk(X)  over 
the  regions  Rj,  specifically: 

11 

Uk(X)  - 2 U^X).  ( A-2) 

i=l 

Then  by  the  symmetry  of  the  regions: 

ul|X>-uk 

l$X)  = u£(X) 

uJ(X)  = Uj(\) 
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uk(x)  = uk°(M 

t 

uJw-uJ'fX) 

and  we  can  write  (A-2)  as 

6 

Uk(X)  = Uj|(X)  + 2 £ U*k(X).  (A-3) 

i=l 

i=£5 

in  order  to  prove  Lemma  3 we  will  need  the  following  relationships  which  are  shown 

in  [1], 

Let 


and 


{bn} 


OO 

n=0 


<cn> 


OO 

n=0 


be  two  sequences  of  complex  numbers  which  are  absolutely  summable  and  let  d be  a com- 
plex number  such  that  0 < |d|  < 1 . We  then  define 


MN  MN  N-n  N-n' 

= I 1 1 ''k'-k‘in-,k'-klck..k+n._n 

1 n=l  n-1  (k,k'€Rj) 

MN  MN  N-n  N-n' 

= S S bk'-kdn  cn-n'  + k'-k. 

N n=l  n'=l  (k,k'€R6) 


( A -4) 


i = 2 5 


( A-5) 


As  N -*•  <»  we  have  the  following  asymptotic  behavior: 


CX1/N) 

r _ . 


q2  = 

N 

” OO 

OO 

r=l 

S=1 

<?3  = 

04  = 


0(1  /N) 
N b0 


Od/N) 

MN 

TTco 


|d|  < I 
d = 1 


|d|  < 1 
d = 1 


ld|  < 1 

d = 1 


( A -6) 


(A-7) 


( A-8) 
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(A-9) 


Qs  = 

% = 


0(  1/N) 
M 

b0C0^ 


N 


0(  1/N) 


M 


_N 

N 


r=l 


ld|<l 

d = 1 

|d|<  1 
d=  1. 


Proof  of  Lemma  2: 

We  have  from  (A-3) 


U,(A)  = U^(A)  + 2 


where 


mN  MN  N-n  N-n' 


U'j(A)  = 


S I S S ^nlC^k+n-^CCtk'+n'-V) 


Wn^lnt,  S'/r, 


• cos  [X(tk+n  - tk)l  cos  lMtk<+n<  - tk')l  • 


(A-10) 


(A-ll) 


In  the  sequel,  U'j(X)  will  be  considered  separately  for  each  region  Rj  and  we  will  use  the  nota- 
tion fn(t)  as  a shortened  version  of  ft  (t). 


Region  R | : 

Since  the  intervals  tk+n  - tk  and  tk%  r.'  - tk'  are  stationary  and  independent  then  we 
can  write: 


UJ(A) = 


1 

(irN/3)^ 


7 7 t 

n 


N-n  N-n' 


( A 1 2) 


where 


-•f 

■'o 


C(t)  cos  (At)  fn(t)  dt. 


Evaluating  the  double  sum  over  (k,k')  in  the  region  Rj  we  have 

MN  MN 

“ “ 


u](A)  = — ~ ^ J znzn'  [(N-n)(N-n'-n+l)-^(N-n)(N-n+l)l(N-n-n'-l)J 


ri*l  n -1 


(A- 13) 
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where  I(  * > is  the  indicator  function,  i.e., 

l(N-n-n'-l)  = 1 if  N-n-n'-I>0 
I(N-n-n'- 1)  = 0,  otherwise. 

As  an  example  ot  evaluating  multiple  sums  over  regions  we  will  go  through  the  steps 
in  evaluating  the  double  sum  in  (A-l  2)  over  R| . 

N-n  N-n'  N-n  N-n' 

XS  1 = S X ,(k<k  + n<k'<k’  + n') 
k.k'eR,  k=l  k'=l 

N-n  N 

" S S Kk'<N-n')l(l  <n)I(n<k'-k)I(l  < n'). 
k=l  k =1 

Now  let  s = k'  - k and  we  have 
N-n  N 

^ ^ l(s  < N - n'  - k)  1(1  < n)  I(n  < s)  1(1  < n') 
k=l  s=l 


N-n-n  N-n'-k  N-n 

= X S 1 = S N-n'-n 


+ 1 - k 


k=l  s=n 


k=l 


= (N-n)(N-n'-n+  I)  - ^ (N  - n)(N  - n + I). 

Now  we  solve  the  integrals  zn  and  zn>  using  Lemma  1 and  (4). 


J0 

Jo 


C(t)  cos  (At)  fn(t)  dt 


e_a|t|  cos(Xt)  0n  —(n  — e^(t'en)  U(t  - en)  dt 


Affn  e^n 
2(n  - 1)! 


e-t[(«+0)-iX]  (t_en)n-l  + c-t[(a+0)+iX]  (t_en)n-l  dt 


■'en 


= Ua  + ft-iM^e^nUa+ft-iX]  r(n) 

+ [(a  + 0)  + iX]~n  e-^n[(a+p)+i\]  r(n) 
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(A-14) 


= Afp  \n  | / e^fla't'iXl  \n~[ 

2 |\l(Of  + 0)  - + \l(a  + 0)  + iA  ] / j 

2 ^(a  + JJ)  - i\|) 


Similarly  we  have  for  zn>: 


.A^"'  /.■*!»-■> I 

n 2 *(((<>+ 01  - iM  j • 

Now  using  (A-14)  and  (A-l  5)  in  (A-13)  we  have 
Mn  Mkj 

I A“  V*  V n+n'  /e-eIa-‘^l  \n  / e-e  I«  — iA  J \n 

U (M  = > > 0n+n  Rtl-H j R nl— .) 

1 4<irN0)2  ^ 4?x  \l(a  + 0)  - iA)/  \l(a  + 0)-i X)/ 


(A-l  5) 


• £(N-n)(N  -n'  -n  + 1)  ~ (N-n)(N  - n + 1 )J  I(N  - n - n'  - 1). 

In  each  ol  the  regions  2-6  we  have  a dependence  between  {C(tk-t-n  - tk) 
cos  A(tk+n  - tk»)  and  {(  < t k'  + n'  - tk')  cos  A(tk'+n'  - tk')}  which  makes  the  expectation  in 
( A-l  1 1 very  ditticult  to  solve.  We  will  deal  with  this  dependency  by  choosing  appropriate 
bounds  for  the  intervals  (tk+n  - tk)  and  (tk'+n'  - tk')  which  will  make  them  independent. 

I tie  bounds  selected  | 1 1 are  as  follows. 

II  7 ^ tjc+fl  - tk,  then  we  have 


•k'-k  ‘=2 

lk+n  “ *k  > = 3 

7>-  tk<-tk  i = 4 

lk+n  " lk  i = 5 

.(tk'-tk)  + (tk+n-tk'+n'>  * = 6 

and  if  6 £ tk'+n>  - tk>  then 

lk'+n'  “ lk+n  ' = 2 

tk'+n"tk+n  » = 3 

h>'  ‘k'+n’  ~ lk'  ' = 4 

lk'+n' " lk'  i = 5 

. lk'+n'  “ lk'  i = 6. 

We  also  will  use  the  bound 
|C(r)  cos  At|  < |C(r)| 
in  evaluating  (A-l  1 ) for  these  regions. 
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Using  these  hounds  we  have  the  following: 


2 1 
U“(M< 


Mn  Mjs|  N_n  N-n'  r°° 


N N N-n  N-n  r°°  r° 

-h-illll  j «t)iV^.)d ,J 

<7r0N)  n=l  n'=l  R,k'€R->  ■' 0 J 0 


|C(t)|  fn'(r)  dr 


( A- 16) 


Uj(M< 


uf(M< 


Uj(M< 


U,(M< 


n -rr 


^ MN  MNN-n  N-n'  /-°° 

(tt/3N)2  n=,  n'=1  k,k'£R  Jq 


{ 


|C(t)|  fn.n'U)  dt  / |C(r)|  fn'(r)  dr 


Region  R-2 

From  (A-l 6),  Lemma  1 and  (5)  we  have 


Ur(A)<—  . 


mn  Mn  N_n  N_n- 


yyyy 

Z,  zl,  zL  (k'-k-i)  / , 

<7rpiN)  n=!  n -1  k,k  €R9  •'e(k'-k) 


e-t(a+0)|t_e(k'-k)|k'- 


/ w 

•-^V-e^n'/  e-t(a+0)  (t  _ en')n'_1  dt 


(n’-l)!-  j , 

•'en 

mn  Mn  N_n  N_n< 

t— » u/K  "K 


I III 


<»0N)-  „=,  n-=|  k k'£R 

dm 

• — ^ e^n  - e-en'(ar+0)  r(  ^ 

(n'-n)!  a + 0 C ""  ’ 


MN  MN  N-n  N-n' 

S I I I wk' 


-k+n'  . 


Now  let  s = k'  - k and  we  have 

2 MN  MN  N-n  N-n'-k 

I Ws+n  l(s  < n - 1 ) I(n  - n'  < s)  I(s  > 1 ). 
,7rpiN)  n=l  n'=l  k=l  s=l 

Consider  the  four  situations: 

( 1 ) n - n'  < 1 and  N-n'-k<n-l 

2 MN  MN  N-n  N-n'-k 

db-S  1 1 1 Ws+n"  [(„-„'<!)  I(N  - n'  - k < n 

U(3N)  n=l  n'=l  k=N-n'-n+l  s=l 

(2)  n - n'  < 1 and  N-n'-k>n-l 
i MN  MN  N-n'-n  n-1 

1111  wStn’  I(n-n'<  l)I(n-  1 < N-n'-k) 


<^N)-  n,-,  k-j 

(3)  n - n'  > 1 and  N-n'-k<n-l 
MN  MN  N-n  N-n'-k 


I I I I — ' I(n -n'>  1)  l(N-n'-k<n 

t v f / . . i 


WN ) ntlnt'lk=N_7-n+li4r-n 

(4)  n - n'  > 1 and  N-n'-k>n-l 


mn  mn  N_n-_n  n_, 

1 1 

’ n=l  n =1  k=l  s=n-n' 


Ws+n  Kn-n'S*  l)I(n-  1 < N-n'-k). 


Summing  expressions  (A- 1 7),  ( A- 1 8),  ( A-J  9),  and  (A-20)  and  simplifying  we  have: 


1 A~ 

1 (7T/JN)2 


IJn  n'  N-n  N-n'-k 

1111 

n-1  l n-1  k=N-n -n+1  s=l 


yys+n 


f I 

n n-1 


M 


N N-n  N-n'-k 


+ ^ 2l  (N-n'-n)Ws+n  + 2^  ^ Yys+n' 

n=n'+l  k=N-n'-n+l  s=n-n' 


n= 1 s=l 


M 


N n-1 


* 1.  1.  (N-n'-n)  Ws+n' 

n=n'+l  s=n-n' 


1 ) (A- 17) 


(A- 18) 


1)  (A-l 9) 


(A-20) 
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The  procedure  for  evaluating  Uj(A)  over  regions  R-3  through  R-6  is  very  similar  to 
its  evaluation  over  R-2.  Therefore,  we  only  state  the  results. 


Region  R-3 


2 MN 

U?(A)<  — y (N  - n')  Wn  (n'  - 1) 

<*0N)~  fa 


Region  R-4 


uj(A)  < 


(7T0N) 


MN  MN  N-n  N-n'-k 

?l  I I I 

.,'=1  „ — I ' 1,  — M ' 


s+n' 


n -1  ln-l+n  k=N-n  -n  s=n-n' 


M 


N n+1 


+ y y (N-n'-n- 1)  Ws+n' 
n=l+n' s=n-n' 


Region  R-5 


M 


N 


uf(X)  < — ^ — - y (N-n)W2n 


(7T0N)“  fa 


Region  R-6 


U*(A)<  — 


, MN 


r-  -> 

~y  (N-n)Wn  ~ + 2-^n 


n=  1 

Finally,  using  (A-3)  and  the  results  for  each  region  we  have: 


, MN 


U,  (A) 


. y (N-n)  W2n 
(TTdN)2^, 


2A2 


Mn  Mn 





1 1 MN 

- (N-n)(N-n+l ) • I(N -n  - n’ - 1 ) + (N  - n')  Wn'(n'  - 1 ) 

n'=  1 
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M 


N 


n'=l 


M 


N N-n  N-n'-k 


M 


N n+1 


2 I I ws+n' + I £ <N-n’-n 


1)  w 


s+n 


n=l+n'  k=N-n'-n  s=n-n' 


n=l+n  s=n-n 


IN  r ->  -i 

^ (N-n)Wn[jy+2-|nJ 


n=l 

Proof  of  Lemma  3 

As  in  the  proof  of  Lemma  2 we  will  use  (A-3)  and  consider  each  region  separately. 

Region  R-l 

We  begin  here  using  (10)  and  noting  that: 


l^(A)<|l^(A)| 


mn  mn  N_n  N_n- 


S 2 S X Etn  |C(tk' ~ lk>  C(lk+n  " V + n 


)l 


n=,  n-=1  kk'€R 


and  using  the  bound  tk'  - tk  > tk+n  - tk,  which  is  independent  of  tk+n  - tk'+n',  we  have: 
MN  MN  N-n  N-n' 

^ ^ ^ ^ E*n  l^^k+n  " *k^Etn  ^(tk+n  ” V+n'^- 
1 „ ' — 1 1 . 1 . 1 


(7r/3N)“  n_,  nT^j  kk'€R 


Substituting  for  fn(t),  the  value  of  C(t)  and  performing  the  integration  results  in: 
MN  MN  N-n  N-n' 

y y y wk'+n-k  whenW=ii^ 

n=l  n =1  k,k  f R , 


'1 


Now  letting  s = k'  - k and  simplifying  we  have 
MN  MN  N-n'-l 


I 

u‘(X)< r 


II  I (N-n'-s)Ws+n'  I(N-n'-n-  1) 


(^N>  n=ln'=l  s=n 
Next  we  consider  regions  R-2  through  R-5  and  beginning  with  ( 10)  we  have: 

MN  MN  N-n  N-n' 


U2(X)=S  111  EtC(tk'-tk)C(tk+n-tk 

n=l  n -1  k,k'€R 


t . t ) 

+n  9 


cos  [X(tk+n  - tk)l  cos  [X(tk-+n-  - tk')l . 


(A-21 
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Using  the  substitutions:  0 = t^,'  - tk,  rj  = tk+n  - t^,  A = - t^'  yields 

Mn  M[sj  j\j_n  N_n'  /•/•/* 

■I  III  C(tf)  C(A)  cos  [ \(6>  + Tj)j  cos  [X(rj  + A)]  (A-22) 

n=l  n'=l  k.k'fRj  JjJ 

' rk'-k<tf  > ‘n-k'+k(T?>  *n'-n+k'-k(A)  dd  dT?  dA- 

Expressing  the  cosines  as  exponentials,  regrouping  terms,  and  using  (1)  we  can  express 
(A-22)  as: 


Mn  Mn  N-n  N-n'  r°° 

— — - y y y y m [ « 

W2H»4l  kTk'€R.  -/o 


e-0(a-iX)  fk,_k(0)  dO  / e-A(a-iX\'_n+k'_k(AklA 


^_k'+k(2X)  + /*  e-A(a+iX)  fn'_n+k'_k(  A)dA  ^-k'+k(0)  . 

Jo 


(A-23) 


This  is  the  same  form  as  (A-4)  which  will  enable  us  to  use  (A-6)  through  (A-9)  and 
hence  determine  the  asymptotic  behavior  of  U‘->(X). 

Region  R-2 

Using  (A-6)  we  can  express  (A-23)  for  region  R-2  as: 


- if 

ln8)  |j=  I -'0  J |_s=l  ■'0 

oo  roo  oo  roo 

y I fjc»_|{(0)  dd  'y  I e_A(o:+iX) 

_r=  I J 0 J |_s=  1 0 


e“A'a-,X)fn._n+k-_k(A)dA 


-A(a+iX)  fn'_n+k'_k(  A>  dA 


u A-24) 


Now,  by  using  Lemma  1 , evaluating  the  integrals,  and  then  simplifying,  we  have  for  ( A-24) 


j)2  N | Z.  \«+(3-iX / 4;  \a+^-iX / 

V / P V c-re(q-iX)  V / A V c-se(a+iX)  I 
4_<  \a+/3-iX/  ZL  \a+/3+iX/ 


.-se(a-iX) 


(A-25) 
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Finally,  summing  ( A-25)  and  rearranging  terms  yields 

U2(A)  = RC  ( f (&-+^~lA) -1 

(tt/3)^  “N  [|j(a  + 0--iA)-0e-e(a_lX)J 

. 21(ort-ff)2  + A2!  - Pc~ea[  2(a+P)  cos  (eX)  - 2X  sin  (eX)|  1 

J(a+/3)2  + A2|  -(3e~ea[2(a+(3)  cos(eX)-2X  sin(eX)]  +02e-2eO! 

as  N -*  °o. 

Equation  (A- 23)  is  solved  in  a very  similar  manner  in  regions  R-3,  R-4,  and  R-5  and 
here  we  will  only  state  the  results. 


Region  R-3 


3 A2 

lAx)  = 

2(  7T/3)2  N 


2 (a+/3)  - 2(k~€01  cos  ( ea ) 


I (a+0)2+A2|  - 2/3e  ea[(a+/3)cos(eX)-Xsin  (eX)|  +(3~e~~ea 


N -*• 

Region  R-4 


A2  MN 


U7(X)=  — RC 

2 2(  7T/3) 2 N 

N -*■  oo 

Region  R-5 

5 A2 

rrm- _ _ do 


CV+/3— iX 

2(ck+/3)  T 

[_(cH-/3-iA)-j3e~e(Q:_ix) 

_(cv+(3 1 2 + A2 

2(n(3)‘ 


(a+0-iA)-1 

2(a+P) 

1 

Ja+P)z  + A^J 

N -+  oo. 


For  region  R-6  we  use  the  substitutions:  0 = tk<  - tk,  r/  = t^  > - tk'  and 
A = lk+n  “ V+n'  ‘n  e4uation  ( A-2 1 ) which  gives  us: 


U2(X)  = 
z 2(tt/3N)2 


MN  MN  N-n  N-n'  rrr 

n=l  n -1  k,k'eRA  JJJ 


C(A)  cos  (X77)  cos  (X(  A + ri  + U )) 


• fk'.iJtf)  fk'+n'-k'(7?>  fk+n-k'-n'(A)  MdAdr). 

Again,  as  in  (A-22)  we  express  the  cosines  as  exponentials  and  regroup  which  yields 
MN  MN  N-n  N-n' 

I I 

C(O)e'X0  fk<_k(0)dtf  (A-26) 


Uj(X)  = 


N "‘N  N-n  N-n'  f r 

'all  11 

(=1  k =1  k.k’tR^  1 J 


2(  7T/3N ) z fa 
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■f 


I'w/c 


C(A)eiXAtk+n_k-_n-(A)dA  M2X)]n'  + J C(tf)e"iXy  fk'_k(0)d0  (A-26) 


■Jc(  A)e^Afk+n_kW(A)dA 

We  note  that  (A-26)  is  of  the  same  form  as  ( A-5),  thus  using  ( A- 1 0)  we  have 


1 MN 


= — = — r« 
2 2(^)2  N 


if 

r=l  J 


+ A lc<V)e~,X°  fT(0)dti  2.  /C(A)e— fs(A)dA. 

s=l 

Using  Lemma  1,  solving  the  integrals  and  simplifying  yields 
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N ->  oo. 


Finally,  using  (A-3)  we  have  the  asymptotic  result 
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N 


(N-n'-s)Ws+n  I(N-n'-n-  1) 


n=l  n'=l  s=n 


[ 
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Proof  of  Lemma  4: 


The  proof  of  this  Lemma  is  similar  to  that  of  Lemma  2.  From  (1 1)  we  have 


U3(X)  = 


1_ 

(tt(3N)2 


MN  MN  N-n  N-n 

ISIS 

,<=i  ,-.'=1  r = i + '=i 


Etn 


^^k'+n'  tk^^k'  *k+n) 


• cos  (X(tk+n  - tk)|  cos  [X(tk'+n>  - tk')J  . 

As  before  we  will  use  (A-3)  and  consider  the  regions  R-l  through  R-6.  Following  the  proce- 
dure used  in  the  proof  of  Lemma  2 we  define  bounds  on  (tk'+n'  - tk)  and  (tk'  - tk+n)  which 
will  make  them  independent  in  each  of  the  six  regions  considered.  We  have  from  ( 1 ],  the 
following  sets  of  bounds. 


For  6 = tk'+n'  ~ *k  we  use  the  b°unds 


(‘k+n-V  + ttk'+n'-V) 

i = 1 

**-k'  - *k*  + *lk'+n'  " lk+n^ 

i = 2 

*lk'+n'  ~ lk+n^ 

i = 3 

V+n'  - V 

i = 4 

lk'+n'  “ lk 

i = 5 

V^k 

i = 6 

and  for  y k tk'  - tk+n  we  simply  use  (tk'  - tk+n)  for  all  six  regions. 

Using  these  bounds  and  (A-3)  we  have  the  following  expressions 
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(A-27) 
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To  solve  the  integrals  and  evaluate  (A-27)  one  follows  the  same  procedure  as  in  the 
proof  of  Lemma  2.  Since  this  procedure  is  quite  lengthy  we  will  state  the  results  for  each  of 
the  six  regions. 
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MN  MN  N-n-n'  N-n'-k 

^ I I "s+"  ■ 
1 1 1 1 


u3<^— !Z<>  2 2 

(TrpN)  0=1,7=!  k=l  s=n 


(A-28) 


Region  R-2 


Mm  M 


N mN  N-n  N-n'-k 


u3W-^-  II  2 I 

(TrpN)-  n=i  n'=n  k=N-n'-n+2  s=l 


MN  MN  n-1  MN  MN  n-1 

+ S S (N-n'-n+1)wS+n,  + ^ ^ ^ (N-n'-n+l ) 1 

n=l  n'=n  s=l  n'=  1 n=n'+l  s=  1 

Mn  Mn  N_n  N-n'-k+l 

+ S S S S WS+"' 

n'=l  n=n'+l  k=N-n'-n+2  s=  1 


( A-29) 


Region  R-3 


U3(X)<— S (N-n')(n'-l)Wn' 
3 (^N)2^, 


(A-30) 


Region  R-4 


U*(A)<-  A".  y Wn  I 2n2N  -3nN  + 3n2-^  - ) 

J (tt/3N)2  1 3 3) 


(A-3 1 ) 


Region  R-5 


2 MN 

U^(X)  < (N-n)W2n 

J (7T0N)2  ^ 


(A-32) 


Region  R-6 


M 


N 


U3(X)  ^ ~ 4 y w"  {^(4  + n(n  - 3))  + n(n(3  - n)  - 2)} 

^ (7T0N)2 


(A-33) 


U3(\)  as 


Now,  using  (A-3)  and  (A-28)  to  (A-33)  and  simplifying  we  can  express  the  bound  for 
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